###########################
### PANDEMICS AND POLITICAL DEVELOPMENT
### THE ELECTORAL LEGACY OF THE BLACK DEATH IN GERMANY
### WORLD POLITICS, VOL. 73, NO. 3 (JULY 2021)
### REPLICATION FILE OF THE EMPIRICAL ANALYSIS
### DANIEL W. GINGERICH & JAN P. VOGLER
###########################

###########################
### REPLICATION FILES PART 3
### PRE-REFORMATION GERMANY
###########################

###########################
### INSTALL REQUIRED PACKAGES
###########################

install.packages("geosphere")



###########################
### DATA MANAGEMENT
###########################

###########################
### DATA ON PARTICIPATIVE ELECTIONS
### BY WAHL (2016)
###########################

data_outcomes=read.csv("wahl_political_data.csv", stringsAsFactors = F)

unique_towns=unique(data_outcomes$city)



###########################
### CREATE NEW VARIABLE: INTRODUCTION OF PARTICIPATIVE ELECTIONS (1300-1500)
###########################

data_outcomes=data_outcomes[which(data_outcomes$century==1300 | data_outcomes$century==1500),]

data_outcomes$ch_elect=NA
for (i in 1:length(data_outcomes$century)){
  if (i/2 == round(i/2))
    data_outcomes$ch_elect[i]=data_outcomes$elect[i]-data_outcomes$elect[i-1]
}

data_outcomes=data_outcomes[which(data_outcomes$century==1500),]



###########################
### DATA ON CONTROL VARIABLES
### BY WAHL (2019)
###########################

data_control_var=read.csv("wahl_2019_data.csv")



###########################
### CREATE FINAL DATASET FOR THE EMPIRICAL ANALYSIS
###########################

data_pre_ref_fin=merge(data_outcomes,data_control_var, by="city", all.x=T)



###########################
### DATA ON BLACK DEATH OUTBREAKS
### BY JEDWAB ET AL. (2019)
###########################

data_bd_outbreaks=read.csv("data_pre_ref_bd_outbreaks.csv", stringsAsFactors=F)

summary(data_bd_outbreaks)



###########################
### CALCULATING BDEI SCORES
###########################

###########################
### COMPUTE THE GLOBAL MAXIMUM DISTANCE
###########################

library(geosphere)
max_store=rep(NA,length(data_pre_ref_fin$century))
for (j in 1:length(data_pre_ref_fin$century)){
  store=rep(NA,length(data_bd_outbreaks$int))
  for (i in 1:length(data_bd_outbreaks$int)){
    store[i]=distm(c(data_pre_ref_fin$long[j], data_pre_ref_fin$lat[j]), c(data_bd_outbreaks$longitude[i], data_bd_outbreaks$latitude[i]), fun = distHaversine)
    max_store[j]=max(store)
  }
}
max_dist=max(max_store)



###########################
### COMPUTE BLACK DEATH EXPOSURE INTENSITY SCORE (V1-5)
###########################

for (j in 1:length(data_pre_ref_fin$century)){
  store=rep(NA,length(data_bd_outbreaks$int))
  for (i in 1:length(data_bd_outbreaks$int)){
    store[i]=distm(c(data_pre_ref_fin$long[j], data_pre_ref_fin$lat[j]), c(data_bd_outbreaks$longitude[i], data_bd_outbreaks$latitude[i]), fun = distHaversine)
    data_pre_ref_fin$BDEI_score[j]=sum(data_bd_outbreaks$int*(1-store/max_dist)^3)
    data_pre_ref_fin$BDEI_score2[j]=sum(data_bd_outbreaks$int*(1-store/max_dist)^6)
    data_pre_ref_fin$BDEI_score3[j]=sum(data_bd_outbreaks$int*(1-store/max_dist)^9)
    data_pre_ref_fin$BDEI_score4[j]=sum(data_bd_outbreaks$int*(1-store/max_dist)^12)
    data_pre_ref_fin$BDEI_score5[j]=sum(data_bd_outbreaks$int*(1-store/max_dist)^15)
  }
}



###########################
### STANDARDIZE BLACK DEATH EXPOSURE INTENSITY SCORE (V1-5)
###########################

data_pre_ref_fin$BDEI_score=scale(data_pre_ref_fin$BDEI_score)[,1]
data_pre_ref_fin$BDEI_score2=scale(data_pre_ref_fin$BDEI_score2)[,1]
data_pre_ref_fin$BDEI_score3=scale(data_pre_ref_fin$BDEI_score3)[,1]
data_pre_ref_fin$BDEI_score4=scale(data_pre_ref_fin$BDEI_score4)[,1]
data_pre_ref_fin$BDEI_score5=scale(data_pre_ref_fin$BDEI_score5)[,1]



###########################
### MAIN EMPIRICAL ANALYSIS IN THE ARTICLE
###########################

###########################
### INTRODUCTION OF PARTICIPATIVE ELECTIONS (1300-1500)
### TABLE 6 IN THE ARTICLE
###########################

lm_ext_01_logit=glm(ch_elect ~ BDEI_score, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_01_logit)

lm_ext_01_logit_bdei2=glm(ch_elect ~ BDEI_score2, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_01_logit_bdei2)

lm_ext_01_logit_bdei3=glm(ch_elect ~ BDEI_score3, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_01_logit_bdei3)

lm_ext_01_logit_bdei4=glm(ch_elect ~ BDEI_score4, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_01_logit_bdei4)

lm_ext_01_logit_bdei5=glm(ch_elect ~ BDEI_score5, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_01_logit_bdei5)



###########################
### INTRODUCTION OF PARTICIPATIVE ELECTIONS (1300-1500) (WITH CONTROLS)
### TABLE 7 IN THE ARTICLE
###########################

lm_ext_ctrl_01_logit=glm(ch_elect ~ BDEI_score + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_ctrl_01_logit)

lm_ext_ctrl_01_logit_bdei2=glm(ch_elect ~ BDEI_score2 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_ctrl_01_logit_bdei2)

lm_ext_ctrl_01_logit_bdei3=glm(ch_elect ~ BDEI_score3 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_ctrl_01_logit_bdei3)

lm_ext_ctrl_01_logit_bdei4=glm(ch_elect ~ BDEI_score4 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_ctrl_01_logit_bdei4)

lm_ext_ctrl_01_logit_bdei5=glm(ch_elect ~ BDEI_score5 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_ext_ctrl_01_logit_bdei5)



###########################
### ADDITIONAL EMPIRICAL ANALYSIS IN THE SUPPLEMENTARY MATERIAL
###########################

###########################
### DESCRIPTIVE SUMMARY STATISTICS
### TABLE A26 IN THE SUPPLEMENTARY MATERIAL
###########################

summary(data_pre_ref_fin)



###########################
### OCCURRENCE OF PARTICIPATIVE ELECTIONS (1500)
### TABLE A27 IN THE SUPPLEMENTARY MATERIAL
###########################

lm_static_ext_01_logit=glm(elect ~ BDEI_score, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_01_logit)

lm_static_ext_01_logit_bdei2=glm(elect ~ BDEI_score2, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_01_logit_bdei2)

lm_static_ext_01_logit_bdei3=glm(elect ~ BDEI_score3, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_01_logit_bdei3)

lm_static_ext_01_logit_bdei4=glm(elect ~ BDEI_score4, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_01_logit_bdei4)

lm_static_ext_01_logit_bdei5=glm(elect ~ BDEI_score5, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_01_logit_bdei5)



###########################
### OCCURRENCE OF PARTICIPATIVE ELECTIONS (1500) (WITH CONTROLS)
### TABLE A28 IN THE SUPPLEMENTARY MATERIAL
###########################

lm_static_ext_ctrl_01_logit=glm(elect ~ BDEI_score + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_ctrl_01_logit)

lm_static_ext_ctrl_01_logit_bdei2=glm(elect ~ BDEI_score2 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_ctrl_01_logit_bdei2)

lm_static_ext_ctrl_01_logit_bdei3=glm(elect ~ BDEI_score3 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_ctrl_01_logit_bdei3)

lm_static_ext_ctrl_01_logit_bdei4=glm(elect ~ BDEI_score4 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_ctrl_01_logit_bdei4)

lm_static_ext_ctrl_01_logit_bdei5=glm(elect ~ BDEI_score5 + elevation + dist_river + romanroad + coast + agsuit + lnpop + ruggednes + urban_potential + trade1300 + protoind1300, family=binomial(link = "logit"), data=data_pre_ref_fin)
summary(lm_static_ext_ctrl_01_logit_bdei5)


